**********************************************************************
*Create scatterplots
**********************************************************************
	*File options 
	local is is 
	local N 500 // N bootstraps

	*Passthrough vars 
	{
	use "$reg/FULLsample`is'", clear
	keep num ye bdtot bdtot_first mcare mcaid paytot2 fte incl_1`is' bought_1`is' wt_cemm ye
	*aha_merge6_1.dta contains hospital-year level data from the AHA survey and 
	*identified mergers, described in Appendix A1, A3. 
	merge 1:1 num ye using "$aha\aha_merge6_1", keepusing(admtot)
	drop if _m==2
	drop _m 
	
	*Prepare hospital-level vars 
	gen     wage = paytot2 / fte 
	gen     comm = 1 - mcare - mcaid
	replace comm = 0 if comm < 0 & !missing(comm)
	drop if wt_cemm==0 
	keep if incl_1`is'
	keep if !bought_1`is' 
	drop paytot2 fte mcare mcaid wt_cemm incl_1`is' bought_1`is'
	
	*Merge to bootstrap weights
	merge 1:1 num_ ye using "$reg\scatter_bweights"
	drop if _m==2
	drop _m 
	
	*Calculate means for each hospital-sample 
	forval i = 1/`N' {
	di "`i'"
	foreach var in bdtot admtot wage comm {
	bysort num: gegen `var'`i' = mean(`var') [fw = bwt_`i']
	}
	}
	foreach var in bdtot admtot wage comm {
	ren `var' `var'0 
	}
	drop ye bwt* 
	
	*Collapse
	bysort num: gen ct = _n 
	keep if ct==1
	drop ct 
	gduplicates report num 
	ren num_prvdr_num _num 
	qui compress 
	tempfile passthru
	save `passthru', replace 
	}	
	
	*Prepare vars 
	use "$sdb\ploteffects", clear		
	{
	*Normalize 
	local _b _b*
	unab _b: `_b'
	foreach var in `_b' {
	qui sum `var', d
	replace `var' = r(p99) if `var' > r(p99) & !missing(`var')
	replace `var' = r(p1)  if `var' < r(p1)  & !missing(`var')
	}
	
	ren _num num_prvdr_num 
	merge m:1 num using "$reg\hets_FULL_`is'"
	drop if _m==2
	drop _m
	ren num_prvdr_num _num 
	
	merge 1:1 _num using `passthru'
	drop if _m==2
	drop _m 

	forval i = 0/`N' {
	qui replace _bcost`i'   = _bcost`i' / 1000
	qui replace _bexptot`i' = _bexptot`i' / 10000
	qui replace _bpaytot`i' = _bpaytot`i' / 10000
	}
	
	*Make room for regression estimates 
	expand 5, gen(new)
	count
	assert r(N) > `N'
	ds 
	foreach v in `r(varlist)' {
	capture noisily replace `v' = .  if new
	capture noisily replace `v' = "" if new
	}
	drop new 
	}
	
	*Plot readmission effect x labor effect 
	*Figure A5(c)
	{
	preserve
	local var _bfte 
	xtile bin_size0 = `var'0, nq(15)
	xtile bin_size1 = `var'0, nq(10)
	xtile bin_size2 = `var'0, nq(5)
	
	local wts 
	local xlab  "Effect on personnel"
	local ylab  "Effect on 90 day readmissions"
	
	*2nd stage bootstrapped regressions 
	reg  _bread0 `var'0 `wts'
	predict pred_read
	local slope = round(_b[`var'`i'],0.001)
	
	sort _num 
	gen B    = . 
	forval i = 1/`N' {
	qui reg  _bread`i' `var'`i' `wts'
	qui replace B = _b[`var'`i'] in `i'
	}
	sum B	
	local se    = round(`r(sd)',0.001)
	
	*Format slope label 
	if abs(`slope') < 1 & `slope' < 0 {
	local slope = abs(`slope')
	local slope = "-0`slope'"
	}
	if abs(`slope') < 1 & `slope' >= 0 {
	local slope = "0`slope'"
	}
	if abs(`se') < 1 {
	local se = "0`se'"
	}
	di "`slope'"
	di "`se'"
	
	*Scatterplot 
	collapse (mean) pred_read `var'0 _bread0 vol, by (bin_size1)
	sum _bread0, d
	local yspot1 = r(p90)
	local yspot2 = r(p25)
	
	twoway (line pred_read `var'0 `wts', lc(black)) (scatter _bread0 `var'0, msym(Oh) mcol(blue*1.5) graphregion(color(white)) xtitle("`xlab'",size(large)) ytitle("`ylab'",size(large)) yla(#5, ang(0)) xla(#6) legend(off) text(`yspot2' 0.8 "Slope: `slope' (`se')", place(e) size(medlarge)) ysize(2.78) xsize(5.5) xlabel(,labsize(large)) ylabel(,labsize(large)))
	graph export "$op\read_`var'.png", as(png) replace
	restore
	}
	
	*Plot expenses effect x acquirer size, NP acquirer size 
	*Figures A3(c,d)
	{
	preserve 
	local var _bexptot
	local het acq_ct 
	*Uncomment for NP only 
	*replace `het' = . if !acq_np 
	
	xtile bin_size0 = `het', nq(15)
	xtile bin_size1 = `het', nq(10)
	xtile bin_size2 = `het', nq(5)
	
	local wts 
	local ylab "Effect on expenses ($10,000/bed)"
	local xlab "Acquirer size (hospitals)"
	
	*Winsorize outliers 
	sum `var'0,d
	replace `var'0 = r(p5)  if `var'0 < r(p5)  & !missing(`var'0)
	replace `var'0 = r(p95) if `var'0 > r(p95) & !missing(`var'0)
	
	*2nd stage bootstrapped regressions
	reg `var'0 `het' `wts'
	predict pred
	local slope = round(_b[`het'],0.001)
	
	sort _num 
	gen B    = . 
	forval i = 1/`N' {
	qui reg `var'`i' `het' `wts'
	qui replace B = _b[`het'] in `i'
	}
	sum B	
	local se    = round(`r(sd)',0.001)
	
	*Format slope label 
	if abs(`slope') < 1 & `slope' < 0 {
	local slope = abs(`slope')
	local slope = "-0`slope'"
	}
	if abs(`slope') < 1 & `slope' >= 0 {
	local slope = "0`slope'"
	}
	if abs(`se') < 1 {
	local se = "0`se'"
	}
	di "`slope'"
	di "`se'"
	
	*Scatterplot 
	collapse (mean) pred `var'0 `het' vol, by (bin_size1)
	sum `var'0, d
	local yspot1 = r(p90)
	local yspot2 = r(p10)
	sum `het', d
	local xspot1 = r(p75)
	drop if missing(`het') | missing(`var'0)
	
	twoway (line pred `het' `wts', lc(black)) (scatter `var'0 `het', msym(Oh) mcol(blue*1.5) graphregion(color(white)) xtitle("`xlab'",size(large)) ytitle("`ylab'",size(large)) yla(#5, ang(0)) xla(#6) legend(off) text(`yspot1' `xspot1' "Slope: `slope' (`se')", place(e) size(medlarge)) ysize(2.78) xsize(5.5) xlabel(,labsize(large)) ylabel(,labsize(large)))
	graph export "$op\het_`var'`het'.png", as(png) replace
	restore
	}
	
	*Plot passthrough 
	*Figure 6(a,b)
	{
	*Calculate total dollar value effects 
	drop _bread* _bpay* 
	forval i = 0/`N' {
	
	*(A) Personnel expense effect 
	*Personnel effect per bed * wage * beds 
	gen a`i' = _bfte`i' * wage`i' * bdtot`i' / 1000
	
	*(B) Commercial revenue effect 
	*Price effect per stay * admits * comm share 
	gen b`i' = _bcost`i' * admtot`i' * comm`i' 
	
	*(C) Operating expenses effect 
	*Expense effect per bed * beds 
	gen c`i' = _bexptot`i' * bdtot`i' * 10 
	drop wage`i' bdtot`i' comm`i' _bexptot`i'
	}
	
	foreach var in a c {
	preserve
	xtile bin_size0 = `var'0, nq(15)
	xtile bin_size1 = `var'0, nq(10)
	xtile bin_size2 = `var'0, nq(5)
	
	if "`var'"=="a" {
	local wts 
	local xlab "Effect on personnel expenses ($'000)"
	local ylab "Effect on commercial revenues ($'000)"
	local loc 7000 12000
	local yscale 
	}
	if "`var'"=="c" {
	local wts 
	local xlab "Effect on operating expenses ($'000)"
	local ylab "Effect on commercial revenues ($'000)"
	local loc 7000 0
	local yscale "yscale(r(-5000 15000))"
	}
	
	*2nd stage bootstrapped regressions 
	reg  b0 `var'0 `wts'
	predict pred_b
	local slope = round(_b[`var'`i'],0.001)
	
	sort _num 
	gen B    = . 
	forval i = 1/`N' {
	qui reg  b`i' `var'`i' `wts'
	qui replace B = _b[`var'`i'] in `i'
	}
	sum B	
	local se    = round(`r(sd)',0.001)
	
	*Format slope label 
	if abs(`slope') < 1 & `slope' < 0 {
	local slope = abs(`slope')
	local slope = "-0`slope'"
	}
	if abs(`slope') < 1 & `slope' >= 0 {
	local slope = "0`slope'"
	}
	if abs(`se') < 1 {
	local se = "0`se'"
	}
	di "`slope'"
	di "`se'"
	
	*Scatterplot 
	collapse (mean) pred_b `var'0 b0 vol, by (bin_size1)
	sum b0, d
	local yspot1 = r(p90)
	local yspot2 = r(p25)
	
	twoway (line pred_b `var'0 `wts', lc(black)) (scatter b0 `var'0, msym(Oh) mcol(blue*1.5) graphregion(color(white)) xtitle("`xlab'",size(vlarge)) ytitle("`ylab'",size(vlarge)) yla(#5, ang(0)) xla(#6) legend(off) text(`loc'  "Slope: `slope' (`se')", place(e) size(vlarge)) ysize(2.78) xsize(5.5) `yscale' xlabel(,labsize(vlarge)) ylabel(,labsize(vlarge)))
	graph export "$op\passthru_`var'.png", as(png) replace
	restore
	}
	}
	